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Abstract. We study an atomic Josephson junction (AJJ) in presence of two 
interacting Bose-Einstein condensates (BECs) confined in a double well trap. We 
assume that bosons of different species interact with each other. The macroscopic 
wave functions of the two components obey to a system of two 3D coupled Gross- 
Pitaevskii equations (CPE). We write the Lagrangian of the system, and from this we 
derive a system of coupled ordinary differential equations (ODE) , for which the coupled 
pendula represent the mechanic analogous. These differential equations control the 
dynamical behavior of the fractional imbalance and of the relative phase of each bosonic 
component. We perform the stability analysis around the points which preserve the 
symmetry and get an analytical formula for the oscillation frequency around the stable 
points. Such a formula could be used as an indirect measure of the inter-species s-wave 
scattering length. We also study the oscillations of each fractional imbalance around 
zero and non zero - the macroscopic quantum self-trapping (MQST) - time averaged 
values. For different values of the inter-species interaction amplitude, we carry out this 
study both by directly solving the two CPE and by solving the corresponding coupled 
pendula equations. We show that, under certain conditions, the predictions of these 
two approaches are in good agreement. Moreover, we calculate the crossover value of 
the inter-species interaction amplitude which signs the onset of MQST. 
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1. Introduction 

The prediction [1] of Bose- Einstein condensation and the reahzation in the laboratory of 
BECs [2] paved the way to a lot of important theoretical and experimental developments. 
One of these is the study of the atomic counterpart [3], HI [5l E] of the Josephson effect 
which occurs in superconductor-oxide-superconductor junctions In [3l HI [5], E], 

realization of AJJ are taken into account from a theoretical point of view. Few years 
ago, Albiez et al. [8] provided an experimental realization of the AJJ. In 2007 Gati et 
al. |9] reviewed the experimental realization of the AJJ focusing on the data obtained 
experimentally with the predictions of a many-body two-mode model [10] and a mean- 
field description. Under certain conditions, a coherent transfer of matter consisting of 
condensate bosons flows across the junction. In the above references the AJJ physics is 
explored in presence of a single bosonic species. The possibility to manage via magnetic 
and optical Feshbach resonances the intra- and inter-species interactions [TTl [T2] makes 
BECs mixtures very promising candidates to successfully investigate quantum coherence 
and nonlinear phenomena such as the existence of self-trapped modes and intrinsically 
localized states. Localized states induced by the nonlinearity were shown to be quite 
generic for multicomponent systems in external trapping potentials. In particular, 
the emergence of coupled bright solitons from the modulational instability of binary 
mixtures of BECs in optical lattices was found numerically in [13]. More sophisticated 
coupled localized states of two-component condensates both in optical lattices and in 
parabolic traps were reported in [H] . The existence of dark-bright states of binary BECs 
mixtures was demonstrated in [15]. On the other hand, the existence of localized states 
of different symmetry type (mixed symmetry states) was numerically and analytically 
demonstrated in [16]. Properties of coupled gap solitons in binary BECs mixtures with 
repulsive interactions were also analyzed in the multidimensional case [17] as well as 
for combined linear and nonlinear optical lattices [H]. Although gap-soliton breathers 
of multicomponent GPE involving periodic oscillations of the two components densities 
localized on adjacent sites of an optical lattice have been found [16] (in analogy to what 
was done for single component case in [6], such states can also be seen as matter wave 
realizations of Josephson junctions), no much numerical and theoretical study has been 
done until now on AJJ of binary mixtures. 

Recently, this has been considered in [191 120] for the case of a bosonic binary 
mixtures trapped in a double well potential, for which a coupled pendula system of 
ODE for the temporal evolution of the relative population and relative phase of each 
component, was derived. Using this reduced system, the authors of [20] have predicted 
the analogous of the macroscopic quantum self-trapping phenomenon for a single bosonic 
component [1]. No comparison between the reduced ODE system and the full GPE 
dynamics has been performed, so that the question of the validity of such prediction 
remains open. For single component condensates, Salasnich et al. [5] have shown that 
a good agreement exists between the results obtained from the GPE and those of the 
ODE. Similar agreement was obtained in [6] for AJJ realized with weakly interacting 
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solitons localized in two adjacent wells of an optical lattice. However, the situation may 
be quite different for multicomponent condensates, due to the interplay of intra- and 
inter-species interactions which enlarge the number of achievable states (for istance, 
mixed symmetry states can exist only in presence of the inter-species interaction) as 
well as their stability, making the system much more complicated. 

The aim of the present paper is just to perform a systematic investigation of possible 
Josephson oscillations which can arise in binary BECs mixtures trapped in a double 
well potential, as a function of the system parameters. To this regard, we derive the 
reduced coupled pendula system proceeding from a Lagrangian formulation and from 
the canonical equations of motion. We show that for certain conditions and range of 
parameters there exists a good agreement between the solutions of the two GPE and 
the predictions provided by the coupled pendula equations. We look for the stationary 
points that preserve the symmetry and study their stability; we get an analytical formula 
for the oscillation frequencies around the equilibrium points. This formula shows the 
possibility to determine the inter-species s-wave scattering length from the frequency. 

We analyze the influence of the inter-species interaction on the temporal evolution 
of each relative population. In particular, by employing the coupled pendula equations 
we show the existence of MQST when the inter-species interaction amplitude is greater 
than a certain value, for which we are able to provide an analytical formula. As done 
by Satjia et al. [20], we calculate the values of the relative populations associated 
to the degenerate GPE states that break the symmetry of the fractional imbalances. 
In addition, we perform the stability analysis by explicitly calculating the associated 
oscillation frequencies. We, moreover, show that the MQST-like evolution obtained by 
solving the coupled pendula equations is close to that one obtained by integrating the 
two coupled GPE. 

Proceeding from the works of Albiez et al. [8] and of Gati et al. [9], we correlate 
our theoretical work with experiments. Finally, we draw our conclusions. 



2. AJJ with two bosonic species: quasi-analytical approach 

We consider two Bose-Einstein condensates of repulsively interacting Bosons with 
different atomic species denoted below by 1 and 2. We suppose that the two BECs 
are confined in a double-well trap produced, for example, by a far off-resonance laser 
barrier that separates each trapped condensate in two parts, L (left) and R (right). We 
assume, moreover, that the two condensates interact with each other. In the mean field 
approximation, the macroscopic wave functions \E'j(r, t), {i = 1,2), of the interacting 
Bose-Einstein condensates in a trapping potential Vtr.ap(r) at zero-temperature satisfy 
the two coupled Gross-Pitaevskii equations 
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Here denotes the 3D Laplacian, and \l/i(r, t) is subject to the normahzation condition 

dxdydz\'^i{x,y,z)\^ = Ni , (2) 




with Ni the number of bosons of the zth species. Similarly, , aj and Qi = Airh'^ai/mi 
denotes the atomic mass, the s-wave scattering length and the intra-species coupling 
constant of the ith species; gij = 2ti1t? aij / rrir {i ^ j) is the inter-species coupling 
constant, with = mimj/{mi + rrij) the reduced mass, and Qij the associated s-wave 
scattering length. In the following, we shall consider both Qi and g^j as free parameters, 



Or 



> 



Q-0.5 




a=1.5; Zjj=2.5 
a=1.2;Zj,=2.5 
a=0.9; Z(,=2.5 



Figure 1. The double well potential (|4]) as a function of z. Lengths are measured in 
units of a^.i = = ax, and energies in units of tnji = huj2 = ^w. 



due to the possibility to change the scattering lengths Oj and aij at will by using the 
technique of Feshbach resonances. Here we take into account the case in which the two 
BECs interact attractively, see [H] and |21j . 

The trapping potential for both components is taken to be the superposition of 
a strong harmonic confinement in the radial (x-y) plane and of a double well (DW) 
potential in the axial (z) direction. For the ith component, we model this trapping 
potential in the form 

,2 



V, 



trap\ 



r) + Vdw{z) 



(3) 



where, for symmetric configurations in the z direction, we model the DW potential as 

2,Z + Zq. 



Vdw{z] 
Vr{z) = -V, 



Vl{z) + Vr{z), Vl{z 

2.Z — Zq 



Sech 



Sech 



a 



a 

Vo = hwAl + Sech^i—)]-^ 

a 



(4) 

i.e. the combination of two Poschl- Teller (PT) potentials, Vl{z) and Vr{z), separated 
by a potential barrier, the height of which can be changed by varying a, centered around 
the points —zq and zq (see Fig. 1). Note that the usage of PT potentials is only for 
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the benefit of improving accuracy in our numerical GPE calculations (see below), taking 
advantage of the integrability of the underlying linear system. We remark, however, that 
the obtained results are of generic validity also for more confining (e.g. not saturating 
to zero at large distances) double well potentials. 

We are interested to study the dynamical oscillations of the populations of each 
condensate between the left (L) and right (R) wells when the the barrier is large enough 
so that the link is weak. To exploit the strong harmonic confinement in the {x-y) plane 
and get the effective ID equations describing the dynamics in the z direction, we write 
the Lagrangian associated to the GPE equations in ([1]) 

^ i=l,2 * 



Vtrar>{r)\^^?--\'^. 



^7..i*ri^.r , (5) 



9l 
2 

where denotes the complex conjugate of and i ^ j and adopt the ansatz [22] 



\ r -|- 

^>i{x,y,z,t) = ^ exp -^r-^ — Mz,t) , (6) 



2a 



where a_Li = \ and the fi(z,t) obey f dz\fi(z)\ = Ni, so that the 

V muJi ~ 
normalization condition Eq. ([2]) is satisfied. By inserting this ansatz ([6]) in ([5]) and 

performing the integration in the radial plane, we obtain the effective ID Lagrangian 

for the fields fi{z,t) 



I dz{\y m^ + ^^)h 

i=l,2 * 

(e.+w(^))i/.p-|i/.r] -~9^,\mfj\') 



(7) 



where the effective parameters for the ID dynamics are given in terms of the original 

= + = 2^' = Trial, + alj ""^^^"^ ^ 

with respect to /j, we obtain the ID GPE for the field ft 
f)f. n'i fff. 

= + + vdw{z) + ^.i/.r + -g^Af^m ■ (8) 

It is possible to study the AJJ dynamics described by Eq. ([8]) by using the two- 
mode approximation discussed by Milburn et al. in [10]. In particular, we assume, for 
each /j, the following time-dependent wave function decomposition 

/,(^,t)=#(t)</>f(^)+^f(t)0f(^), (9) 

where 

= v^e^'?W , (10) 
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with a=L,R, and a constant total number of particles given by N^" + = I'^/'f + 
\ijf^{t)\^ = Ni, with J^^ dz\(f)2{z)\'^ = 1 and J^^ dz(j)f{z)(j)f{z) = 0. Neglecting terms 
of order greater than two in the overlaps of the 0's, we can write the Lagrangian in 
terms of A^" and as 



L= ^^i^^ - ^^^^^ - ^^^^ - ^^^^ 



=1,2 



iVM^cosf 



R\2\ 



where 



K, -- 



dz 



2mj dz 
dz (6?)\ U- 



+ {Vdw + 



2m,- ai 



12 



~g,2 / dzm 



a\2/ Aa\2 
2) 1 



— I dz 



-2mi dz dz 

(12) 

One may get a good approximation for the functions (f)f{z) and 4>^{z) when the 
double well potential Vdw{z) is such that the two lowest energy eigenvalues of the 
corresponding Schrodinger equation constitute a closely spaced doublet well separated 
from the higher excited levels, and the ^'s are not too large (see, for example, [10]). If 



the real symmetric function (t)f{z) and the real antisymmetric function 0f (2) are the 
wave functions of the ground state and of the first excited state, respectively, then (^f (2;) 
and (t)f{z) may be chosen as 

4>!{z) - <Pt{z) 



(13) 



Remember that 



and 



dz\<pf{z)\'' = I ^nd dz<Pi 



V2 

satisfy the relations: dz\(j)f{z)\'^ = 

[z) = 0. Having chosen Vdw{z) as the sum 
of two of two Poschl- Teller (PT) wells (see Eq. (jlj)), the functions (f)^{z) and 4'f'{z) 
may be analytically calculated following a perturbative approach. Let us consider the 
eigenvalues problem corresponding to Eq. ([8]) with cji = cjij = 0. We know exactly the 
wave functions for this eigenvalues problem when the potential is given by a single Va{z) 
[a = L,R), for example Vl{z). The wave function of the ground state is [23] 



AL,PT), 



A[l - Tanh'^ 



.z + Zq 



B,/2 



1 
2 



2mjVoa^ 1 



(14) 



In Eq. (fT4|) A, equal for both sides, ensures the normalization of the wave function. 
Since we are assuming that the two lowest energetic levels are well separated from the 
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higher ones, when the potential is perturbed by the presence of Vji{z), we look for 
the eigenstates in the form of a linear superposition of 4>l^'^^^ and For each 

component, to the first order of such a perturbative theory, the ground state wave 
function ch^iz) and the first excited state wave function (hf^iz) read 



Ma{ 



AL,PT)i 



f '(^)) 



(15) 

Here Ms and Ma ensure the normalization of (pf{z) and (pf{z), respectively. Note that 
Ki is equal to {Ef — E^)/2, with and Ef, also perturbatively calculated, the energies 
associated to (pf{z) and (pf{z), respectively. The quantity 



()f (z) and (Pfiz) 
(Ef - Ef) 



00 



h 



(16) 



is the Rabi frequency. This frequency characterizes the oscillations of a particle between 
the states and (j)f-. By using Eq. (fT5l) in the decomposition (fT3l) . we are able to write 
the functions 4>f{z) and (f)f{z) in terms of (f)[^'^'^\z) and (j)\^'^'^\z) in the following way 



{Ms + Ma) 



,{L,PT) 



{Ms - Ma)<P. 



{R,PT) 



V2 

{Ms - M^)0f '^^^ + {Ms + M^)#'^^) 
^ V2 



{17) 



Note that 4>f{z) and ^^(-z), and the associated energies, may be numerically found as 
the wave functions of the two lowest states of the eigenvalues problem corresponding 
to Eq. ([S]) in absence of interactions. Then, by using the decomposition fll3p . one cal- 
culates the functions 4>f{z) and (t)f{z). We have verified that the perturbative theory 
provides practically the same results as the numerical approach . 

Let us, now, focus on the Lagrangian (fTTl) . The conjugate moments of the 
generalized coordinates and hd^ are given by 

p^. = ^ = , p,. = = -iVf . (18) 

The Hamiltonian of the system is 

^ = - E [P^f + M^^"^ - E ^K,y^iw^cos{ef - ) + 

i=l,2 4=1,2 
U L jjR 

+ E [y + i-^ef ] + U^2Pe^^Pe^ + U^2Pe^pgH . (19) 

1=1,2 

The evolution equations for the fractional imbalance Zi = {N^ — N^) / Ni and for the 
relative phase 9i = Of — Of for each component are derived from the canonical equations 
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associated to the Hamiltonian 0191) 

By subtracting the equation for from the one for pgi^, we obtain the equation for 
the temporal evolution of Zi. Similar arguments leads to the equation for 6i(t). In the 
following, we shall assume two wells to be symmetric, i.e. -Ef = E^, = = Ui, 
= = Ui2. The fractional imbalance and the relative phase for each component 
vary in time according to the following (coupled pendula) equations: 



9i = — — + cos9i + . 

^ hy^l - zf ^ 

(21) 

Note that when U12 = 0, equations (12T|) reduce to the usual equations for a single 
component obtained in [1]. At this point, we observe that it is possible to obtain Eq. 
(I2T]) proceeding from the following equations 



where H is 



H=-Y1 '^KiN.i^il - zf) cos^,] + Y^^' + U,2NiN2Z^Z2. (23) 

i=l,2 2=1,2 

To fix ideas, let us consider, as done in [12] , a mixture of two bosonic isotopes of the same 
atom, so that one may have different intra-species and inter-species interactions. For 
simplicity, for the time being we will neglect the mass difference between the two species. 
We compare the temporal behavior of the fractional imbalances obtained integrating Eq. 
(IH]) with the Ziit) obtained by solving the coupled differential equations ( 12T|) . To this 
end, in the wave functions ([Hj) we assume that the quantities Bi are the same for both 
species, and in the double well potential (jlj) we set a = 1.2 and zq = 2.5, as illustrated by 
the continuous line of Fig. 1. We calculate the parameters Ui, U12, and Ki = K2 = K 
by using Eqs. ( fT2l) . Note that the value of K (equal, for the above values of a and 

zq, to 0.0148), provided by the last formula of Eqs. (fT2l) coincides with — . Here tq 

To 

is the oscillation period of Zi{t) obtained by numerically solving the ID GPE ([8]) when 
Ui = U12 = 0. This period is just equal to 27t/uj with u given by Eq. fITB]) . In Fig. 2, 
the first panel of each Zi{t) graph shows the perfect agreement between the two coupled 
ID GPE and the coupled pendula equations when the bosons do not interact at all. In 
the second panel of each Zi{t) of Fig. 2, U is finite and U12 = 0. Finally, in the last 
panel, Ui and Uu are both finite. From the second and the third panels of Fig. 2, we see 
that the solutions of the ODE system ( l2Ti) with K = 0.0148 (dot-dashed lines) shows 
a certain displacement with respect to the ones (dashed lines) predicted by solving the 
two ID GPE (IH]). The continuous lines, in the second and in the third panels of Fig. 2, 



U 
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represent the solutions of the ODE system fl21l) with K obtained via a fitting procedure. 
The best-fit K^s that we found are equal to 0.0155 and 0.0151 when only f/j 7^ and 
when Ui and Uu are both finite, respectively. Notice that these values are very close to 
the above theoretical estimate of K, so the coupled pendula equations may be used to 
consistently describe the AJJ physics. 

To gain physical insight in the behavior of the system, we carry out the stability 
analysis. We study the oscillations around the points for which the time derivatives of 
Zi and of 6i are zero, i.e. around the stationary points. We diagonalize the Jacobian 
matrix associated to Eqs. (l2Tll . The eigenvalues A allows us to determine the frequencies 
of oscillation around the points of equilibrium. If the A's are of the form iu (with u a 
pure real number), u will be the oscillation frequency around stable equilibrium points. 
By performing the above analysis, we obtain the following classes of stationary points: 

(i) z, = ,Z2 = ,61 = ,62 = 

(ii) ^1 = , ^2 = , ^1 = 0, ^2 = Tvr 

(iii) Zi = 0,Z2 = 0, 9i = TTT, 02 = 

(iv) zi = 0,Z2 = 0, 9i = =F7r, 62 = =Fvr, 
zi = 0, 22 = 0, 61 = =F7r, 62 = ±7T . 

We observe that the points of a given class are characterized by the same A. The 
stationary solution (i) represents a stable equilibrium if Uu satisfies the condition 

- A(i) < f/12 < A^') 
^(1) ^ Vi2K, + N^U,)i2K2 + N2U2) 



provided 



2K 

U.>-^. (25) 



The small amplitude oscillations frequency, say uj^^\ around the point (i) is 
= l{Ki{2K\ + UiNi) + K2{2K2 + U2N2) ± A)'/' , 



AK^K2{4.K^K2 + 2K1U2N2 + 2K2UiN^ - UI2N1N, 



2 

+ U1U2N1N2) + {Ki{2Ki + UiNi) + ^2(27^2 + U2N2)fj (26) 

with + and — corresponding to the normal modes of the linearized system associated to 
Eq. ([HI). When 21 (0) = ^^2(0) < 1, f/i = f/2, i^i = i^2, and A^i = A'a, the fractional 
imbalances Zi oscillate around the point (i) according to the law 

Zi{t) = Zi{0) cos u^^k . (27) 

Let us operate in Eqs. (^^, (125|1 . and (126|1 the replacement Ui — f/j. Then, we 
obtain the stability conditions and the oscillation frequency associated to the class (iv). 
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Figure 2. Fractional imbalance of the two bosonic species vs. time. Here A^i = 100 
and N2 = 150. In each plot of Ziif), from top to bottom: U\ ^ Ui = U12 = 0; 
Ui = U2 ^ 0mi,Ui2 = 0; C/i = C/2 = C/ = 0.001, ;7i2 = -U/2. The dashed line 
represents data from integration of Eq. ([8]), the dot-dashed line represents solution 
of Eq. (|2ip with Ki = K2 = K = 0.0148, and the continuous hne, in the second 
and third panels, represents solution of Eq. pT|) with the best- fit A''s, say K^f. 
For each z,{t), Ktf = 0.0155 {Ui = U2 = 0.001, C/12 = 0) and Ki,f ^ 0.0151 
{Ui ^ U2 = U ^ 0.001, C/12 = ~U/2). We used the initial conditions zi(0) = 0.1, 
2:2(0) — 0.15, and 0^(0) = 0. Time is measured in units of = (^2)"^ = ^ 

lengths are measured in units of a^.i = a^^2 ^ and energies in units of Hlj. 



We have verified that for the stationary points of type (ii), under certain conditions 
(analytically achievable but very complicated), the eigenvalues of Jacobian matrix are 
all of the form iu. The small amplitude oscillations frequency, say u;''^-', around the point 
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11 IS 



J^^ = ^{M2K^ + f/iiVi) + K2{2K2 - U2N2) ± Af^ , 

- AKiK2iAKiK2 - 2K1U2N2 + 2K2U1N1 + UI2N1N2 
UiU2N^N2) + (Ki(2Ki + f/iiVi) + K2{2K2 - U2N2)f ^'^ 



f28) 



For the oscillations of Zi around the point (ii), one may use arguments analogous to the 
ones employed for the class (i). If we start, now, from the points of the class (ii), and 
replace Ui with — t/j, we get the stability regions and the oscillation frequency for the 
point of type (iii). Let us focus, to fix the ideas, on the frequency fl26l) . and on the 
formula of Eq. f[T21) which gives the inter-species interaction amplitude U12. We note 
that gi2 is directly related to the inter-species s-wave scattering length. This quantity, 
then, can be determined from the oscillation frequency (1261) once one keeps fixed Ki, Ui, 
and Ni. We will discuss this point with more details in the following. 

At this point, it is worth observing that, because of the non linearity associated to 
the inter- and intra-species interactions, there is a class of degenerate GPE eigenstates 
that breaks the Zi symmetry. Let assume that Ui = U2 = U, Ki = K2 = K, and 
Ni = N2 = N. In correspondence of 6i = ir, we have looked for non zero stationary 
solutions of the system (12T1) . We have found four classes of fractional imbalances 
corresponding to the Zi broken symmetry; we have verified that two of these classes 
do not correspond to a stable equilibrium. Let us consider the two classes describing 
stable equilibrium, say / and //. For the class /, we have 

SI) - 




2K 

.N{U + Ui2) 

^2,s6 — Alb 1 (29) 

provided \ {U + Ux2)\ > 2K/N. When < f/ < 2K/N, the solution ^ is always stable, 
and the corresponding oscillation frequency is 



{I) 



^\j{N{U + U12)) -4:KK (30) 



For U > 2K/N, the solution ([29]) is stable when 
Ui2>uS = ^(^-UN/2K 



31/3 _ {gUN/2K + v/3 + 81(f/A^/2K)2)2/3 
32/3(9f/iV/2K + v/3 + 81(t/A^/2A')2)V3 

(31) 

The corresponding oscillation frequency is 

(I cr) 

It is possible determine the crossover value, say U{2 , of the inter-species 
interaction strength signing the onset of the self-trapping. We start evaluating the 
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Hamiltonian ( !23l) in Zi = z^^]f^ and 6i = vr. Let us denote the energy obtained in this way 
by E^^\ This energy reads 

E^'^ = -—— + N\U + U^,) . (33) 

U + Ui2 

Then, we evaluate the Hamiltonian (123!) at t = 0, i.e. H{zi{0),9i{0)). We require that 

H{z,{0),e,m>4KN (34) 

with 4KN the value got by E'^^^ when z^^]fj = 0, i.e. when 2K = N{U + Uu). Then, we 
get 



{I,cr) 



K 



UN 



iV^i (0)^2(0) 2K ^^^^^ 
+ 2 E v/l-^i(0)2cos^,(0)] . (35) 

i=l,2 

(I cr) 

When the condition Uu > U12 is satisfied, the system will be self-trapped. For the 
solution of the type //, we have 



h,sb - ±1 



2K 

^N{U-Ui2)- 

43 = > (36) 

provided that \{U-Ui2)\ > 2K/N. When <U < 2K/N, this solution is always stable 
and is characterized by the oscillation frequency 



ill) 



^\j[N{U -U,2)) -AK^ . (37) 



When U > 2K/N, the solution ([MD is stable if 
Uu<uii'^ = ^{uN/2K 



31/3 _ {-9UN/2K + y/3 + 81{UN/2K)'^ 



\2/3 



3-^/3 {-9UN/2K + v/3 + 81(f/iV/2if)2)i/3 
and the corresponding oscillation frequency reads 



(3^ 



-■^■-y{m-u4.4Ki^). (39, 

Also for the solution //, it is possible to determine the inter-species interaction 
amplitude which signs the self-trapping onset, say u[2^'^^\ by using the same argument 
employed for the class /. Also in this case, we proceed by evaluating the Hamiltonian 
( !23|) in Zi = zl^2 — ^- Let us denote the energy obtained in this way by E^^^\ 

This energy reads 
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Figure 3. Fractional imbalances of the first and second bosonic species as functions vs. 
time. Here Ni = 100 and N2 = 150, Ki = A'2 = 0.0148, and Ui = U2 = U ^ O.lKi. 
Plots are for different values of Uu- U12 = —U/20 and U12 = —U/2 (the two upper 
panels of each Zi{t), from left to right); U12 = —U and U12 = —1.2U (the two lower 
panels of each Zi{t), from left to right). We used the initial conditions zi(0) =0.1, 
^2(0) = 0.15, and 9i{0) = 0. Here, we report the profiles and I/2P as functions of 
z and for different values of time, as displayed in the figure. The quantities are 



obtained by integrating Eq. 
2. 



when U12 = —1.2U. Units are as in Fig. 1 and Fig. 
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Figure 4. Fractional imbalance of the first and second bosonic vs. time. The dashed 
line represents data from integration of Eq. ([8]), the dot-dashed line represents solution 
of Eq. pT|) with Ki = K2 = K = 0.0148, and the continuous line represents solution of 
Eq. jUl) with the best-fit K equal to 0.0151. Here is Ui = U2 = O.IK, U12 = -O.UK, 
Ni ^ N2 = 100. We used the initial conditions Zj(0) = 0.4 and 9,(0) = 0. The units 
are as in Fig. 1 and in Fig. 2. 



The condition to find f/jg '^'^ ^^e Eq. f l34|) . and we obtain that t/jg ''^^ coincides with 

An interesting task, now, is to analyze the influence of the interaction between the 
two BECs on the temporal evolution of the fractional imbalances. We have flxed V\ 
and f/2, and analyzed z\{f) and 2^2 (^) for different values of the inter-species interaction 
amplitude U\2- The greater is the absolute value of f/12, the greater is the deformation 
of the oscillations around (^(t)) = 0, as shown in the two upper panels (from left to 
right) and in the flrst lower panel (from the left) of each Zi(f) represented in Fig. 3. 
Note that as long as the oscillations are harmonic (see the flrst panel of each ^i(t)), the 
time evolution of the fractional imbalances may be described in terms of a carrier wave 
of frequency uJc-, given by Eq. fl26l) with U\2 = 0, modulated by a wave of frequency 
ijjra- The frequency ujc is much greater than ujra- We found that there exists a value of 
the inter-species interaction amplitude for which the relative population in each trap 
oscillates around a non zero time averaged value, {zit)) ^ 0, which corresponds to the 
macroscopic quantum self-trapping (MQST) as discussed in [1] for a single component. 
To support this interpretation, we have studied the behavior of the density proflles of 
the two species as function of z and for different values of time. In particular, to flnd 
|/j(z,t)P, we have numerically solved the two coupled GPE ([8]) for those values of the 
interaction amplitudes for which the self-trapped is predicted to occur by the coupled 
pendula equations. We have summarized the results of this analysis in the last two pan- 
els of Fig. 3. Finally, it is interesting to observe that there is a good agreement between 
MQST predicted by the coupled pendula equations ( 12T|) and the one got by numerically 
solving the two ID GPE ([8]). This comparison is displayed in Fig. 4, where we can 
observe oscillation around a non zero time averaged value of fractional imbalances. We 
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followed the same fitting procedure adopted in obtaining Fig. 2. 

At this point, we observe that the question of the AJJ with two bosonic species could 
be addressed also from the experimental point of view. This could be done considering, 
for example, BECs binary mixtures of two bosonic isotopes of the same alkali atom 
|12j . The works referenced in [8] and [9] provide the ideal guide lines for this kind of 
experiments. By engineering a double well potential as suggested, for istance, by Gati 
et al. [9], the measured fractional imbalances Zi{t) could be compared with the ones 
obtained both by solving the Eq. ([8]) (see [8] and [9]) and the ODE system ( |2T]) . It 
could be possible to measure, moreover, the inter-species s-wave scattering length ai2 
by using, to fix the ideas, the frequency (l26l) . The mixture could be prepared in such 
a way that both all the conditions ^i(O) = 22(0), Ki = K2, Ui = U2, Ni = N2 - see 
the discussion about the four classes (i)-(iv) of stationary points - are verified and to 
have small amplitude oscillations around the stationary point (i). Each Zi oscillates 
according to the law 2;j(0) cos(u;i^'*t) (see Eq. ( |2711 ). Let us suppose to fix both the 
characteristic quantities (in our case, they are uji, a, Zq) of the trapping potential - the 
group of Heidelberg displayed how this is possible for a given class of double well traps 
[8], [9] - and the intra-species s-wave scattering length Oj. Then, the functions 0° are 
known. Then, the Eqs. (fT2!) provide the parameters Ki and f/°. The measure of the 
period of Zi{0) cos{uj^yt) leads to the corresponding frequency 00^^, i.e. the left-hand 
side of Eq. fl26l) . The solution of this equation gives the parameter Ui2- By the mean 
of the second line of Eq. ( fT2l) one gets gu, and, then a^. 

3. Conclusions 

We have analyzed the atomic Josephson effect in presence of a binary mixture of BECs. 
We have written the Lagrangian of the system, from which we have derived a system 
of coupled differential equations which governs the dynamical behavior of the fractional 
imbalance and of the relative phase of each component. We have analyzed the stable 
points that preserve the symmetry, and we have obtained an analytical formula for the 
frequency oscillations around these equlibrium points. To this regard, one of the most 
interesting features is the possibility to know the inter-species s-wave scattering length 
from these frequencies. We have shown that in correspondence of precise values of the 
inter-species interaction amplitude, the relative populations oscillate around a non zero 
time averaged value. This behavior corresponds to MQST, a well-known phenomenon 
when only one component is taken into account. We have compared the predictions of 
CPE with the ones of the coupled pendula equations. We have performed this compari- 
son in the case of total absence of interaction, in the case in which only the intra-species 
interaction is present, and in the case in which also the inter-species interaction is in- 
volved. We have found that, under certain conditions, the predictions of CPE agree 
with those ones of the coupled pendula equations. We have shown that, under certain 
hypothesis, it is possible to obtain analytical expressions for the inter-species interac- 
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tion amplitudes which signs the onset of the self-trapping. Finally, we have commented 
about the possibility to correlate our theoretical work with the experiments proceeding 
from the works of the group of Heidelberg, see [8] and [9] . 
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